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Abstract 


A  numerical  simulation  has  been  performed  of  turbulent  flow  near  and  through  a  two- 
dimensional  array  of  rectangular  buildings  (bluff  bodies)  in  a  neutrally  stratified  back¬ 
ground  flow.  The  model  used  for  the  simulation  was  the  steady-state  Reynolds-averaged 
Navier-Stokes  equations  with  an  isotropic,  linear  eddy  viscosity  formulation  for  the  Reynolds 
stresses;  and,  the  eddy  viscosity  was  determined  using  a  high-Reynolds  number  form  of  the 
k-e  turbulence-closure  model  with  the  boundary  conditions  at  the  wall  obtained  with  the 
wall-function  approach.  The  resulting  system  of  partial  differential  equations  along  with  the 
concomitant  boundary  conditions  for  the  problem  was  solved  using  the  SIMPLE  algorithm 
in  conjunction  with  a  non-orthogonal,  colocated,  cell-centered,  finite  volume  procedure. 
The  predictive  capabilities  of  the  numerical  simulation  of  urban  flow  are  validated  against  a 
very  detailed  and  comprehensive  wind  tunnel  data  set.  Vertical  profiles  of  the  mean  stream- 
wise  and  vertical  velocities  and  the  turbulence  kinetic  energy  are  presented  and  compared 
to  those  measured  in  the  wind  tunnel  simulation.  In  addition  to  the  comparison  between 
modeled  and  measured  results,  the  numerical  simulation  is  used  to  aid  in  the  interpretation 
of  certain  features  of  the  flow  observed  in  the  wind  tunnel  simulation.  It  is  found  that  the 
overall  performance  of  the  model  is  good — most  of  the  qualitative  features  in  the  disturbed 
turbulent  flow  field  in  and  near  the  building  array  are  correctly  reproduced,  and  the  quan¬ 
titative  agreement  is  also  generally  fairly  good  (especially  for  the  mean  velocity  field).  The 
underestimation  of  the  turbulence  energy  levels  may  arise  from  the  inability  of  the  model 
to  capture  the  subtle  effects  of  secondary  strain,  especially  those  associated  with  curvature, 
on  the  turbulence. 
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Executive  summary 


Introduction:  It  is  anticipated  that  Canadian  Forces  (CF)  in  the  foreseeable  future  will 
have  to  fight  in  or  protect  urban  areas,  whether  in  battle,  peace-making,  peace-keeping, 
or  counter-terrorist  operations.  The  increased  awareness  and  importance  accorded  by  the 
public  worldwide  and  their  governments  to  maintain  appropriate  defences  against  chemical 
and  biological  warfare  (CBW)  agents  in  an  urban  (built-up)  environment,  the  prediction 
of  casualties  and  human  performance  degradation  resulting  from  such  releases,  and  the 
development  of  operational  procedures  and  regulations  to  control,  mitigate,  and  monitor  the 
fate  of  CBW  agents  in  urban  areas  with  high  population  densities,  will  require  mathematical 
modeling  of  urban  wind  flows  and  dispersion.  In  this  regard,  it  should  be  noted  that  the 
prediction  of  urban  canopy  flow  is  in  principle  pre-requisite  to  or  co-requisite  with  the 
prediction  of  scalar  (e.g.,  CBW  agent)  dispersion  within  a  canopy. 

Results:  A  numerical  simulation  has  been  performed  of  turbulent  flow  near  and  through  a 
two-dimensional  array  of  rectangular  buildings  (bluff  bodies)  in  a  neutrally  stratified  back¬ 
ground  flow.  The  model  used  for  the  simulation  was  the  steady-state  Reynolds-averaged 
Navier-Stokes  equations  with  an  isotropic,  linear  eddy  viscosity  formulation  for  the  Reynolds 
stresses;  and,  the  eddy  viscosity  was  determined  using  a  high-Reynolds  number  form  of  the 
fc-6  turbulence-closure  model  with  the  boundary  conditions  at  the  wall  obtained  with  the 
wall-function  approach.  The  resulting  system  of  partial  differential  equations  along  with  the 
concomitant  boundary  conditions  for  the  problem  was  solved  using  the  SIMPLE  algorithm 
in  conjunction  with  a  non-orthogonal,  colocated,  cell-centered,  finite  volume  procedure. 
The  predictive  capabilities  of  the  numerical  simulation  of  urban  flow  are  validated  against  a 
very  detailed  and  comprehensive  wind  tunnel  data  set.  Vertical  profiles  of  the  mean  stream- 
wise  and  vertical  velocities  and  the  turbulence  kinetic  energy  are  presented  and  compared 
to  those  measured  in  the  wind  tunnel  simulation.  In  addition  to  the  comparison  between 
modeled  and  measured  results,  the  numerical  simulation  is  used  to  aid  in  the  interpretation 
of  certain  features  of  the  flow  observed  in  the  wind  tunnel  simulation.  It  is  found  that  the 
overall  performance  of  the  model  is  good — most  of  the  qualitative  features  in  the  disturbed 
turbulent  flow  field  in  and  near  the  building  array  are  correctly  reproduced,  and  the  quan¬ 
titative  agreement  is  also  generally  fairly  good  (especially  for  the  mean  velocity  field).  The 
underestimation  of  the  turbulence  energy  levels  may  arise  from  the  inability  of  the  model 
to  capture  the  subtle  effects  of  secondary  strain,  especially  those  associated  with  curvature, 
on  the  turbulence. 

Significance  and  Future  Plans:  The  fc-e  turbulence-closure  model  with  an  isotropic, 
linear  eddy  viscosity  is  the  simplest  complete  turbulence  model  (viz.,  no  advance  knowledge 
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of  any  property  of  the  turbulence  is  required  for  the  simulation  other  than  the  initial  and/or 
boundary  conditions  for  the  problem)  that  is  available  currently,  and  its  moderately  good 
predictive  performance  of  the  developing  flow  through  an  array  of  two-dimensional  build¬ 
ings  without  the  need  to  adjust  any  closure  constants  is  encouraging.  This  model  may  be 
useful  as  a  general-purpose  simulator  of  urban  flows  since  it  is  simple  enough  to  be  tractable 
numerically  and,  hence,  not  require  excessive  computing  time.  However,  before  this  model 
(or,  variants  of  it)  can  be  used  for  this  purpose,  it  needs  to  be  validated  against  more  ex¬ 
perimental  data  for  flow  in  various  obstacle  arrays.  From  this  perspective,  the  availability 
of  extensive  and  accurate  experimental  data  for  flow  through  various  obstacle  arrays  is  ur¬ 
gently  required.  Finally,  the  accurate  prediction  of  urban  dispersion  requires  knowledge  of 
the  mean  wind  and  turbulence  through  the  urban  area.  The  utility  of  the  urban  flow  simu¬ 
lation  models  investigated  here  for  provision  of  the  disturbed  wind  field  statistics  required 
for  a  physically-based  dispersion  model  needs  further  investigation. 


Yee,  E.,  Lien,  F.  S.,  and  Cheng,  Y.  2001.  Numerical  Simulation  of  Flow  Near  and  Through  a 
Two-Dimensional  Array  of  Buildings.  DRES  TR  2001-115,  Defence  Research  Establishment 
Suffield. 
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Introduction 


In  recent  years,  a  great  effort  has  been  made  to  understand  the  impact  of  large  numbers  of 
discrete  bluff  obstacles  (e.g.,  buildings,  trees,  and  other  obstructions  within  a  small  domain 
inside  a  large  city)  on  the  atmospheric  boundary-layer  flow  because  of  its  importance  in 
many  aspects  of  meteorology,  wind  engineering,  and  environmental  science.  In  particular, 
prediction  of  how  the  speed,  direction,  and  turbulence  of  the  wind  changes  as  it  blows 
over  urban  surfaces  is  in  principle  pre-requisite  to  the  prediction  and  understanding  of  the 
transport  and  turbulent  diffusion  of  contaminants  within  an  urban  environment.  Indeed, 
in  response  to  increasing  concerns  on  the  use  of  chemical  and  biological  warfare  (CBW) 
agents  against  military  and  civilian  populations  in  urban  areas,  there  has  been  a  concerted 
effort  by  a  number  government  and  military  research  laboratories  and  agencies  to  develop 
models  that  can  predict  the  turbulent  flow  and  dispersion  in  urban  areas  (e.g.,  Defense 
Threat  Reduction  Agency  Urban  Dispersion  Model  Working  Group,  Department  of  Energy 
Chemical  and  Biological  National  Security  Program). 

As  far  as  numerical  modeling  of  urban  flows  is  concerned,  most  of  the  effort  has  been 
concentrated  on  modeling  the  flow  around  isolated  structures  or  single  buildings  rather 
than  around  groups  of  obstacles  or  building  clusters.  Prognostic  models  that  solve  the 
Reynolds- aver  aged  Navier- Stokes  equations  have  been  used  to  describe  the  flow  field  around 
an  isolated  cuboid  or  building  (Patterson  and  Apelt  [1];  Zhang  et  al.  [2];  Zhou  et  al.  [3]; 
Cowan  et  al.  [4])  or  the  interaction  of  two  buildings  (Lee  [5]).  These  models  used  the 
simple  k-e  turbulence  model  in  the  Reynolds-averaged  Navier-Stokes  (RANS)  equations  to 
determine  the  mean  flow  and  turbulence  and,  hence,  were  basically  nonlinear  flow  models 
solved  using  either  finite  difference,  finite  element,  or  finite  volume  solution  techniques. 
Computational  fluid  dynamics  (CFD)  models  have  also  been  applied  to  describe  the  flow 
within  a  specific  urban  structure  (e.g.,  urban  street  canyon).  In  particular,  Hunter  et  al. 
[6]  have  used  the  fc-e  model  of  turbulence  to  describe  the  three-dimensional  characteristics 
of  flow  regimes  within  an  urban  canyon. 

Generally  speaking,  urban  flow  and  dispersion  within  various  arrays  of  obstacles  are  ad¬ 
dressed  on  a  case-by-case  basis.  From  a  theoretical  point  of  view,  the  problem  appears  to 
be  complex  and  unwieldy  and  there  appears  to  be  no  idealized  starting  point  for  analysis. 
From  an  applications  point  of  view,  it  is  easy  to  get  lost  in  the  complexities  of  particular 
cases  without  ever  having  to  see  that  there  are  common  structures  underlying  the  entire 
subject.  However,  from  a  numerical  viewpoint,  it  appears  that  the  numerical  simulation 
models  based  on  RANS  equations  are  general  enough  to  be  applied  mutatis  mutandi  to  the 
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problem  of  disturbed  flow  within  building  clusters  where  mutual  flow  interferences  between 
buildings  need  to  be  taken  into  account.  The  principal  problem  with  the  application  of 
such  CFD  models  to  arrays  of  buildings  is  that  there  is  currently  very  little  benchmark 
mean  flow  and  turbulence  data  that  could  be  used  to  compare  with  the  model  predictions. 
In  this  regard,  it  should  be  stressed  that  CFD  modeling  applied  to  a  particular  problem 
always  yields  an  answer.  For  example,  we  could  apply  a  simplistic  turbulence  model  to 
a  complicated  problem  and  get  an  answer.  It  may  be  a  poor  answer,  but  it  is  an  answer 
nonetheless,  and  there  is  no  “alarm”  automatically  warning  the  user  that  a  particular  turbu¬ 
lence  closure  may  be  applied  outside  its  realm  of  applicability.  For  this  reason,  it  is  crucial 
that  before  any  CFD  model  with  its  attendant  turbulence  closure  is  applied  to  a  particular 
class  of  problems  (e.g.,  flow  through  an  array  of  buildings),  that  the  model  predictions  are 
subjected  first  to  critical  and  comparative  judgements  against  high-quality  data  sets  before 
such  model  predictions  are  taken  as  “ground  truth”. 

In  this  report,  we  use  the  standard  k-e  model  of  turbulence  with  a  linear  eddy  viscosity 
model  to  simulate  the  separated  flow  associated  with  a  two-dimensional  array  of  surface- 
mounted  obstacles  set  within  a  turbulent  boundary-layer  under  neutrally-stratified  flow.  We 
show  a  detailed  comparison  with  a  wind  tunnel  experiment  of  the  average  flow  fields  and 
turbulence  kinetic  energy  for  this  obstacle  array.  These  results  are  of  more  than  theoretical 
interest  since  they  provide  information  that  can  be  used  to  interpret  the  data  from  this 
wind  tunnel  experiment  and  to  ascertain  basic  phenomenological  aspects  of  the  experiment 
which  are  not  apparent  from  the  laboratory  data. 
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The  wind  tunnel  experiment 


The  wind  tunnel  experiment  is  fully  described  by  Brown  et  al.  [7]  and  only  the  important 
details  of  the  experiment  will  be  presented  here.  The  experiments  were  conducted  in  an 
open-return  type  wind  tunnel  at  the  U.  S.  Environmental  Protection  Agency’s  Fluid  Model¬ 
ing  Facility.  This  wind  tunnel  has  a  working  test  section  of  length  18.3  m,  width  3.7  m,  and 
an  adjustable  roof  of  height  approximately  2.1  m  to  eliminate  streamwise  pressure  gradients 
and  allow  for  a  non- accelerating  free  stream  flow. 

A  schematic  of  the  two-dimensional  building  array  used  for  these  wind  tunnel  experiments 
is  shown  in  Figure  1.  This  array  consisted  of  seven  rectangular  blocks,  each  block  with  equal 
height,  H,  and  length,  L,  of  150  mm.  The  blocks  were  spaced  H  apart  in  the  streamwise 
direction.  The  streamwise,  cross-stream,  and  vertical  velocity  components  were  measured 
with  a  pulsed  wire  anemometer  using  a  pulsing  rate  of  10  Hz  and  a  sampling  time  of  120  s 
at  each  measurement  location.  The  pulsed  wire  anemometer,  unlike  a  standard  X- array 
hot-wire  anemometer,  can  be  used  to  measure  velocities  in  high  intensity  turbulence  with 
reversing  flows  such  as  that  which  exists  within  an  array  of  obstacles. 

Measurements  of  the  three  components  of  the  instantaneous  velocity  were  made  at  1016 
coordinate  positions.  These  measurements  consisted  of  81  vertical  profiles  of  mean  velocity, 
turbulence  velocity  variance,  and  turbulence  kinetic  energy.  These  profiles  extended  from 
3.5H  upstream  of  the  first  building  to  7.5H  downstream  of  the  last  building.  Each  vertical 
profile  consisted  of  either  8  (over  the  rooftops)  or  16  (over  the  urban  canyons,  and  either 
upstream  or  downstream  of  the  array)  measurement  positions  up  to  a  vertical  height  of  3H 
above  the  ground  surface.  All  streamwise  distances  are  measured  relative  to  the  location  of 
the  upstream  wall  of  the  first  building  at  x/i7  =  0. 
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Numerical  simulation 


Governing  equations,  turbulence  closure,  and  boundary  conditions 


For  flow  near  and  through  an  array  of  obstacles,  we  solve  the  Reynolds-averaged  equations 
for  continuity  and  momentum  for  steady,  incompressible  flow.  These  equations  are  as 
follows: 

P:  =  o,  (1) 

(yOC<i 


dUi  ,, 


d(u,Uj)  ap  e  f  au,  \ 

Here,  Ui  and  u'  are  the  mean  and  fluctuating  velocities  in  the  Xj-direction,  respectively,  P 
is  the  kinematic  pressure,  is  the  kinematic  fluid  viscosity,  and  (■)  denotes  the  ensemble¬ 
averaging  operation.  In  Equations  (1)  and  (2),  the  convenient  Einstein  summation  con¬ 
vention  is  assumed,  which  prescribes  that  if  any  of  the  indices  is  repeated  in  a  term,  a 
summation  over  that  index  is  implied.  Subscripts  1,  2,  and  3  refer  to  directions  along  x 
(streamwise)  ,  y  (spanwise),  and  ^  (vertical),  respectively;  averaged  quantities  are  in  capital 
letters  or  are  designated  by  angle  brackets.  The  presence  of  the  Reynolds  stresses, 
in  the  mean  momentum  equations  [Equation  (2)]  implies  that  the  latter  are  not  closed, 
viz.  they  contain  more  unknown  quantities  than  there  are  equations.  Closure  requires  that 
some  approximations  be  made  in  prescribing  the  Reynolds  stresses  in  terms  of  the  mean 
flow  quantities.  The  most  popular  statistical  turbulence  closures  are  linear  eddy  viscos¬ 
ity/gradient  transport  models  that  assume  a  linear  relationship  between  turbulent  stresses 
and  mean  velocity  gradients: 


,  ,  ,,  fdUi  dUj\  2,^ 


where  ut  is  the  kinematic  eddy  (turbulent)  viscosity  and  k  =  is  the  turbulence 

kinetic  energy  (TKE).  Substituting  Equation  (3)  into  Equation  (2)  leads  to  the  following 
closed  form  for  the  momentum  equation: 

d{UiU^)  on  d  /  ^dUA  , 


where 


U  =  P+-k 
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is  the  augmented  pressure. 


The  kinematic  eddy  viscosity,  Ut,  is  on  dimensional  grounds  the  product  of  a  turbulence 
length  scale.  A,  and  a  turbulence  velocity  scale,  uq.  In  the  standard  k-e  turbulence  model, 
the  turbulence  velocity  and  length  scales  are  determined  as  uq  =  Cll'^k'^l'^  and  A  =  u^jt 
where  e  is  the  dissipation  of  turbulence  kinetic  energy  and  is  a  closure  constant.  Com¬ 
bining  the  turbulence  length  and  velocity  scales,  we  obtain  the  familiar  expression  for  the 
eddy  viscosity 


vt  =  C^,  —  . 


The  standard  model  uses  the  following  transport  equations  for  k  and  e: 


d{lJ^k)  _  d 


i/t\  fdUi  dUAdUi 

^  C7fc/  dxj_  dxi  )  dxj 


dxj 


de]  ^  e  fdUi  dUAdUi  ^ 


The  equations  contain  five  closure  constants;  namely,  C^,  Uh,  Cu,  and  €2^.  The  standard 
k-e  model  employs  values  for  constants  that  are  arrived  at  by  comprehensive  data  fitting 
over  a  wide  range  of  turbulent  flows  [8]: 


=  0.09;  (Tfc  =  1.00;  a,  =  1.30;  Cu  =  1.44;  C2,  =  1.92.  (7) 

The  turbulent  diffusion  of  k  and  e  in  Equations  (6a)  and  (6b),  respectively,  are  represented 
using  a  gradient  diffusion  hypothesis  with  the  Prandtl  numbers  cr*,.  and  used  to  connect 
the  eddy  diffusivities  of  k  and  e  to  the  eddy  viscosity  Vf 


Since  the  computational  domain  is  chosen  to  be  large  compared  with  the  two-dimensional 
array  of  obstacles,  the  flow  at  its  boundary  is  not  affected  by  any  of  the  obstacles.  Hence, 
free  boundary  conditions  are  imposed  at  all  air-to-air  boundaries  in  the  flow  domain.  At 
the  inflow  boundary,  the  measured  distributions  of  streamwise  mean  velocity,  U  =  Ui,  the 
vertical  mean  velocity,  W  =  1/3  =  0,  and  the  turbulence  kinetic  energy,  k,  are  used.  This 
information  is  obtained  from  the  vertical  profiles  of  U  and  k  measured  furthest  upstream 
from  the  upstream  face  of  the  first  building  (viz.,  at  x/H  =  -3.5,  or  x  =  -500  mm,  where 
X  =  0  corresponds  to  the  upstream  wall  of  the  first  building).  Above  z  =  3H,  measurements 
for  U  and  k  are  not  available;  and,  for  the  inflow  boundary  conditions,  the  vertical  profile 
of  U  above  2  =  3H  was  assumed  to  adhere  to  a  power-law  profile  with  an  exponent  of 
0.16  (corresponding  to  a  nominal  neutral  boundary  layer  simulated  in  the  wind  tunnel 
using  spires  [7]),  and  the  vertical  profile  of  k  was  assumed  to  be  constant  above  z  =  3H 
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(corresponding,  as  such,  to  a  deep  constant  stress  layer).  Unfortunately,  at  the  inflow,  e 
was  not  measured  and,  consequently,  estimated  simply  as  e  =  where  the  dissipation 

length  scale  I  was  chosen  as  approximately  1.2H,  The  length  scale  I  is  unknown  at  the 
inflow  boundary,  but  variation  of  its  value  from  O.IH  to  1.2H  showed  that  it  had  only  a 
very  small  impact  on  the  final  solution.  Far  downstream,  at  the  outflow  boundary,  we  set 
dUfdx  =  dk I dx  =  defdx  =  W  =  0.  At  the  upper  boundary,  we  set  dU  jdx  —  W  =  0, 
dkfdz  =  defdz  =  0. 


At  all  solid  boundaries  (ground,  obstacle  walls,  obstacle  roofs)  in  the  computational  domain, 
wall  functions  are  used  for  mean  velocities  and  turbulence  quantities.  The  implementation 
of  the  wall  boundary  conditions  starts  with  the  evaluation  of 


n 


(8) 


where  np  is  the  normal  distance  of  the  near  wall  node  P  to  the  solid  surface  (ground, 
wall,  or  roof),  is  the  shear  stress  at  the  wall,  and  p  is  the  fluid  density.  It  is  assumed 
that  the  first  grid  point  closest  to  the  solid  surface  is  located  in  the  logarithmic  region, 
30  <  np,  where  the  flow  is  in  local  equilibrium  (viz.,  the  production  is  exactly  balanced  by 
the  dissipation).  In  this  case,  =  \/t^/p  is  given  by 


and 

=  iln(£;n+),  (96) 

where  k,  is  von  Kar man’s  constant  {k  ^  0.41  for  smooth  and  rough  surfaces),  Vt^p  is  the  mean 
velocity  parallel  to  the  wall  at  the  near  wall  node  P,  np  =  u^npfp  is  the  dimensionless 
distance  from  the  wall  at  the  near  wall  node  P,  and  E  is  the  wall  roughness  parameter 
which  for  a  hydrodynamically  smooth  wall  assumes  the  value  E  ~  9.793  [9]  (obtained  from 
correlation  of  measurements  for  smooth  surfaces).  In  the  current  simulations,  all  the  solid 
surfaces  are  assumed  to  be  smooth.  The  wall  boundary  condition  for  the  tangential  velocity 
requires  a  relationship  between  the  velocity  at  the  first  grid  point  above  the  wall  and  the 
wall  shear  stress,  which  in  accordance  to  Equations  (9a)  and  (9b),  can  be  expressed  as 


Finally,  the  TKE  dissipation  rate  at  the  first  grid  point  P  above  the  wall  can  be  set  to  the 
value 
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Numerical  method 


Our  computational  domain  for  simulating  the  two-dimensional  array  of  buildings  used  in  the 
wind  tunnel  experiment  (cf.  Figure  1)  extended  in  the  streamwise  direction  from  x/H  —  —5 
io  x/H  =  15,  with  the  upstream  wall  of  the  first  building  at  x  =  0.  The  height  of  the 
domain  was  QH.  The  results  of  the  numerical  simulation  reported  below  were  obtained 
with  an  orthogonal  251  x  51  grid  shown  in  Figure  2.  The  grid  lines  were  concentrated  near 
the  solid  surfaces  (ground,  rooftop,  and  building  walls)  and  the  spacing  between  grid  lines 
was  gently  stretched  at  increasing  distances  from  the  solid  surfaces. 

The  numerical  method  employs  a  structured,  non-orthogonal  (allowing  boundary-fitted  co¬ 
ordinates),  fully  colocated  (viz.,  all  flow  variables  are  stored  at  the  same  location),  cell- 
centered,  finite  volume  approach  for  the  discretization  of  the  computational  domain.  Within 
this  non-orthogonal  system,  the  velocity  vector  is  decomposed  into  its  Cartesian  compo¬ 
nents,  and  these  are  the  components  to  which  the  momentum  equations  relate.  Advective 
volume  face  fluxes  are  approximated  using  the  upwind  discretization  scheme,  and  physical 
diff'usive  cell  fluxes  are  approximated  using  a  conventional  second-order  central  differencing 
scheme.  The  SIMPLE  algorithm  [10]  was  used  for  pressure  correction.  Here,  mass  conti¬ 
nuity  is  enforced  by  solving  a  pressure  correction  equation  which,  as  part  of  the  iterative 
sequence,  steers  the  pressure  toward  a  state  at  which  all  mass  residuals  in  the  cells  are 
negligibly  small.  In  conjunction  with  the  colocated  grid  use  here,  this  method  is  known  to 
provoke  chequerboard  oscillations  in  the  pressure  field,  reflecting  a  state  of  pressure- velocity 
decoupling.  To  avoid  this,  the  widely  used  method  of  Rhie  and  Chow  [11]  is  adopted  to 
non-linearly  interpolate  the  cell  face  velocities  from  the  nodal  values  (at  the  center  of  the 
cells).  The  interpolation  scheme  essentially  introduces  a  fourth-order  pressure  diffusion 
correction  that  smooths  out  the  pressure  if  it  oscillates  rapidly.  In  the  numerical  scheme, 
the  transport  equations  for  the  mean  velocity,  turbulence  kinetic  energy,  and  viscous  dis¬ 
sipation  rate  and  the  pressure-correction  equation  are  solved  sequentially  and  iterated  to 
convergence,  defined  by  reference  to  the  sum  of  the  absolute  cell  residuals  for  mass  and 
momentum  components.  Stone’s  strongly  implicit  procedure  [12]  was  applied  to  solve  the 
system  of  algebraic  equations  arising  from  the  discretization  of  the  transport  equations  and 
pressure-correction  equation.  Iterative  refinement  of  all  fields  was  continued  until  the  sum 
of  the  absolute  cell  residuals  for  mass  and  for  momentum  was  satisfied  to  within  0.1%  of 
the  total  mass  and  momentum  fluxes,  respectively,  at  the  inflow  boundary.  For  the  current 
simulation,  this  required  about  780  outer  iterations  which  took  approximately  210  s  on  an 
800  MHz  Pentium  III  PC. 
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Results  and  discussion 


Mean  streamwise  velocity 


Figure  3  compares  model  predicted  vertical  profiles  of  the  mean  streamwise  velocity  U  at 
six  fixed  x-locations  covering  the  range  —2.5  <  x/H  <  —0.1  in  the  impingement  zone 
directly  upstream  of  the  first  building.  At  x  =  -225  mm  and  -150  mm,  the  wind  speed 
reduction  below  about  z  =  2H  relative  to  the  far  upstream  profile  is  well  predicted  by 
the  model,  although  the  wind  speeds  below  about  zjH  =  1/3  are  slightly  over-estimated. 
At  X  =  —75  mm  and  —38  mm,  the  wind  speed  reduction  below  and  wind  speed  increase 
above  the  level  z/H  4/3  relative  to  the  upstream  profile  is  reproduced  well.  However, 
the  reverse  flow  near  the  surface  at  these  x-positions  is  underestimated  by  the  model.  In 
particular,  streamfunction  contours  of  the  recirculation  zone  (not  shown)  in  front  of  the  first 
building  indicate  that  the  predicted  reattachment  length  of  this  zone  occurs  at  z/iJ  —0.33, 
whereas  the  measurements  suggest  that  the  reattachment  occurs  at  about  z/H  ^  —0.52. 
Close  to  the  upstream  wall  of  the  first  building  at  x  ==  -15  mm,  the  significant  wind  speed 
reduction  below  the  building  height  H  owing  to  the  adverse  pressure  rise  (horizontal  pressure 
gradient)  upstream  of  the  solid  obstacle  is  modeled  well,  as  is  the  moderate  speed-up  above 
the  building  top. 

Figure  4  shows  vertical  profiles  of  the  mean  streamwise  velocity  U  at  five  fixed  x-locations 
(0.1  <  x/H  <  0.9)  over  the  roof  of  the  first  building.  Generally  speaking,  the  speed-up  over 
the  roof  compared  with  the  reference  upwind  profile  and  the  presence  of  a  thin  separation 
zone  near  the  roof  surface  with  flow  reversal  is  well  matched  between  the  predictions  and 
the  measurements,  although  the  magnitude  of  the  reversed  flow  is  underestimated  slightly. 
Figure  5  shows  reasonable  agreement  between  the  measurements  and  the  model  predictions 
of  U  at  various  fixed  x-locations  within  the  first  urban  canyon  (1.1  <  x/H  <  1.9).  An 
important  feature  captured  by  the  model  predictions  is  the  very  strong  shear  layer  that 
forms  immediately  downstream  of  the  rear  face  of  the  first  building,  whose  signature  is 
revealed  in  the  inflection  point  in  U{z)  at  or  near  the  building  height  H.  Note  that  the 
high  values  of  mean  shear  dU/dz  just  above  building  height  over  the  first  urban  canyon 
(Figure  5)  are  reproduced  well  by  the  model.  The  speed-up  above  about  z/H  ^  4/3  is 
predicted  correctly  over  the  first  urban  canyon,  although  the  magnitude  of  the  reverse  flow 
in  the  canyon  below  about  z/H  ^  2/3  is  underestimated  slightly.  The  reverse  velocity  on 
the  underside  of  the  standing  vortex  within  the  canyon  is  comparable  to  the  forward  velocity 
at  the  canyon  top,  and  this  feature  appears  to  be  reproduced  by  the  model  predictions. 
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Figure  6  compares  computed  and  measured  mean  streamwise  velocity  U  at  five  x-stations 
over  the  rooftop  of  the  second  building  (2.1  <  x//f  <  2.9).  Again,  wind  reduction  below  and 
moderate  speed-up  above  about  z/H  ~  4/3  is  well  predicted.  The  model  predicts  that  no 
separation  occurs  over  the  second  rooftop  (in  contrast  to  the  behavior  over  the  first  rooftop) 
which  accords  well  with  the  observations.  Vertical  profiles  of  U  at  five  x-locations  within 
the  second  urban  canyon  (3.1  <  x/H  <  3.9)  and  over  the  rooftop  of  the  third  building 
(4.1  <  x/H  <  4.9)  are  displayed  in  Figures  7  and  8,  respectively.  Note  that  the  predicted 
magnitude  of  the  wind  shear  over  the  second  urban  canyon  is  larger  than  what  is  observed 
in  the  wind  tunnel  measurements.  The  primary  reason  for  this  discrepancy  appears  to  be 
the  reduced  vertical  turbulent  diffusion  of  shear  stress  in  the  high  shear  layer  bordering  the 
urban  canyon  top — a  defect  manifesting  itself  by  the  greater  level  of  shear  strain  in  this  area 
than  observed.  The  reduced  vertical  turbulent  diffusion  results  from  the  underestimation  of 
the  turbulence  kinetic  energy,  which  implies  an  underestimation  in  the  turbulent  viscosity. 
Because  of  this,  the  high  shear  layer  generated  near  the  urban  canyon  top  dissipates  more 
slowly  in  the  streamwise  direction  along  the  canyon  and  the  succeeding  rooftop  than  in  the 
wind  tunnel  (cf.  Figures  7  and  8).  In  consequence,  the  prediction  of  speed-up  above  and 
wind  reduction  below  the  level  z/H  ^  4/3  is  slightly  over-predicted  and  under-predicted, 
respectively. 

Vertical  profiles  of  both  the  computed  and  measured  mean  streamwise  velocity  U  within 
the  third  urban  canyon  and  the  fourth  rooftop  are  almost  identical  to  the  profiles  obtained, 
respectively,  for  the  second  urban  canyon  (Figure  7)  and  the  third  rooftop  (Figure  8).  In 
consequence,  U  appears  to  have  reached  streamwise  equilibrium  (i.e.,  the  mean  streamwise 
velocity  field  is  fully  developed)  by  the  second  canyon,  a  feature  that  both  the  wind  tunnel 
measurements  and  the  numerical  simulation  agree  on.  Within  the  array  interior  (x  >  3H), 
it  appears  that  the  mean  streamwise  velocity  within  the  second  urban  canyon  and  above 
third  rooftop  should  be  “typical”  of  their  neighbors  (viz.,  there  is  little  variation  in  the 
amplitude  of  the  streamwise  wind  speed  modulation  from  one  building  to  the  next),  and 
that  in  the  array  interior  a  periodic  boundary-condition  might  have  been  used  to  model  the 
flow  in  a  single  representative  building  and  urban  canyon. 

Figures  9  and  10  compare  computed  and  measured  mean  streamwise  velocity  U  in  the  exit 
region  of  the  array  directly  downstream  of  the  last  building  covering  the  range  13.1  <  x/H< 
20.5.  The  vertical  profile  of  17  at  x  =  1965  mm  lies  just  15  mm  downstream  of  the  rear  wall 
of  the  last  building.  Below  the  building  height  H,  the  vertical  profile  of  U  at  this  x-location 
exhibits  essentially  no  shear  (i.e.,  dU /dz  ^  0)  and  this  absence  of  shear  is  reproduced  well 
by  the  model  predictions  as  is  the  high  value  of  mean  shear  at  or  near  the  building  height 
H.  The  model  predicts  the  existence  of  a  recirculation  zone  behind  the  last  building  with  a 
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reattachment  point  at  x  ~  2475  mm,  which  is  in  good  conformance  with  the  observations. 
Further  downstream,  in  the  far  wake  region,  the  vertical  profiles  of  U  are  recovering  towards 
the  far  upwind  reference  state.  It  is  seen  that  the  model  underpredicts  slightly  the  rate  of 
momemetum  recovery  in  this  region.  It  is  interesting  to  note  that  even  at  the  furthest 
downstream  measurement  location  (at  xfH~  20.5  or  at  x/H  =  7.5  downstream  of  the  last 
building),  the  streamwise  mean  velocity  has  still  not  recovered  to  its  upwind  equilibrium 
state. 


Mean  vertical  velocity 


Comparisons  of  measured  and  modeled  vertical  profiles  of  the  mean  vertical  velocity  W , 
at  six  fixed  x-positions  covering  the  range  —2.5  <  x/H  <  —0.1  in  the  impingement  zone 
directly  upstream  of  the  first  building,  are  shown  in  Figure  11.  As  the  front  face  of  the 
first  building  is  approached,  the  vertical  velocity  is  positive  as  air  is  forced  up  and  over  the 
building.  This  feature  in  the  observations  is  reproduced  well  by  the  modeled  results,  as  is 
the  position  and  magnitude  of  the  maximum  positive  vertical  velocity  at  or  near  the  top  of 
the  building  (viz.,  at  z/H  1).  For  x  >  -38  mm,  the  downward  (negative)  vertical  velocity 
below  about  z/H  ^  2/3,  resulting  from  the  flow  hitting  the  front  face  of  the  obstacle  and 
flowing  downwards  to  create  the  region  of  reversed  flow  in  front  of  the  first  building,  is 
underestimated,  especially  closer  to  the  front  face  (e.g.,  compare  measured  and  computed 
results  for  x  =  -15  mm  in  Figure  11).  Indeed,  in  front  of  the  first  building,  the  model 
predicts  that  there  is  a  vertical  pressure  gradient  with  a  high  pressure  region  centered  at 
about  z/H  ^2/3  that  forces  the  air  over  the  first  building  above  this  level  and  down 
towards  the  ground  below  this  level.  This  prediction  appears  to  be  in  good  conformance 
with  the  measurements  of  the  mean  vertical  velocity. 

Profiles  of  measured  and  computed  mean  vertical  velocity  W  covering  the  range  0.1  < 
x/H  <  0.9  over  the  rooftop  of  the  first  building  are  given  in  Figure  12.  The  upward  motion 
of  air  above  the  first  building  rooftop  agrees  well  with  the  wind  tunnel  measurements.  The 
thin  separation  zone  and  negative  vertical  velocity  in  this  zone  is  captured  in  the  numerical 
simulation.  It  should  be  noted  that  the  streamline  curvature  over  the  first  building  rooftop 
as  well  as  immediately  upstream  of  this  rooftop  is  concave  upwards,  implying  a  destabilising 
(exhange-enhancing)  influence  on  the  turbulence.  This  feature  of  the  flow  curvature  is 
reproduced  in  the  model  predictions.  Figure  13  compares  computed  and  observed  vertical 
profiles  of  W  for  five  fixed  x-stations  (1.1  <  x/H  <  1.9)  within  the  first  street  canyon. 
The  prediction  that  the  downward  (negative)  vertical  velocity  on  the  downstream  end  of 
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the  street  canyon  is  larger  in  magnitude  than  the  upward  (positive)  vertical  velocity  on 
the  upstream  end  of  the  canyon  is  in  conformance  with  the  observations.  However,  the 
magnitude  of  both  the  upward  and  dov/nward  vertical  velocities  are  under-predicted  by  the 
model,  with  the  degree  of  under-prediction  being  larger  closer  to  the  canyon  walls. 

Measured  and  modeled  vertical  profiles  of  W  over  the  second  building  rooftop  are  exhibited 
in  Figure  14  for  x-stations  in  the  range  2.1  <  x/H  <  2.9.  In  contrast  to  the  positive 
vertical  velocity  over  the  first  building  rooftop,  there  is  a  weak  downward  vertical  motion 
over  the  second  rooftop  arising  from  the  slight  downward  curvature  (implying  a  stabilizing 
or  exchange-suppressing  influence)  of  the  mean  flow  streamlines  above  the  upstream  edge  of 
the  second  building.  This  feature  in  the  flow  is  well  reproduced  in  the  numerical  simulations. 
Figures  15  and  16  display,  respectively,  vertical  profiles  of  W  over  the  second  street  canyon 
(3.1  <  x/H  <  3.9)  and  the  third  building  rooftop  (4.1  <  x/H  <  4.9)  in  comparison  with 
the  simulation  using  the  present  model.  With  reference  to  Figures  13  and  15,  it  is  seen 
that  strength  of  the  vortex  in  the  second  urban  canyon  is  smaller  than  that  in  the  first 
canyon  in  the  sense  that  vertical  motions  are  stronger  in  the  first  canyon  than  in  the  second 
canyon.  This  feature  is  captured  in  the  numerical  simulation,  although  the  magnitude  of 
the  upward  motion  on  the  upstream  side  and  downward  motion  on  the  downstream  side 
is  under-predicted  in  the  second  canyon.  However,  the  degree  of  under-prediction  is  not 
as  severe  as  it  was  for  the  first  street  canyon  (cf.  Figure  13).  Nevertheless,  the  locations 
of  the  peak  upward  and  downward  vertical  velocities  in  the  street  canyons  appear  to  be 
well  predicted.  The  mean  flow  streamlines  over  the  third  building  rooftop  are  essentially 
horizontal  and  there  is  negligible  vertical  mean  velocity  over  this  rooftop.  This  aspect  of 
the  flow  is  represented  well  in  the  numerical  simulations  (cf.  Figure  16). 

Vertical  profiles  of  both  the  computed  and  measured  mean  vertical  velocity  W  within  the 
third  street  canyon  and  over  the  fourth  building  rooftop  are  identical  (approximately  or 
better)  to  the  profiles  obtained,  respectively,  within  the  second  street  canyon  (Figure  15) 
and  the  third  rooftop  (Figure  16).  Hence,  in  conformance  with  the  behavior  of  U,  the  mean 
vertical  velocity  field  appears  to  be  fully  established  by  the  second  urban  canyon  within 
the  array  interior,  a  feature  of  the  flow  that  is  in  good  overall  conformance  with  the  model 
predictions. 

Figures  17  and  18  compare  computed  and  measured  profiles  of  mean  vertical  velocity  W 
in  the  exit  region  of  the  array  directly  downstream  of  the  last  building  covering  the  range 
13.1  <  x/H  <  20.5.  Generally  speaking,  computed  and  measured  W  over  the  range  of  the 
exit  region  for  which  observations  are  available  differ  by  less  than  about  5%,  except  within 
the  recirculation  bubble  at  x  =  1965  and  1988  mm.  In  particular,  the  upward  vertical 
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motion  just  downstream  of  the  last  building  in  the  recirculation  zone  is  overestimated  by 
about  80%.  The  prediction  of  the  weak  downward  vertical  motion  just  upstream  of  the 
reattachment  point  is  in  good  conformance  with  the  measurements.  Finally,  the  rate  of 
recovery  of  the  vertical  mean  momentum  beyond  the  reattachment  point  in  the  far  wake 
region  to  its  far  upstream  reference  state  is  reproduced  accurately. 


Turbulence  kinetic  energy 


Figure  19  presents  vertical  profiles  of  turbulence  kinetic  energy,  A:,  in  the  impact  region 
upstream  of  the  first  building  at  the  same  x-stations  as  the  mean  velocity  field  described 
above  (—2.5  <  xjE  <  —0.1).  Although  the  increase  in  k  above  the  building  height  H  as 
the  building  is  approached  is  reasonably  well  predicted,  the  positions  of  the  maximum  k  for 
^  ^  (i.e.,  in  the  region  of  recirculation  zone  just  ahead  of  the  front  face  of  the  first 

building)  are  too  far  away  from  the  ground  surface.  This  feature  seems  to  be  accompanied 
by  a  second  defect  in  the  model  predictions  (perhaps  indirectly  provoked  by  the  first  in  the 
upward  shift  of  the  position  of  maximum  turbulence  generation) — namely,  an  insufficient 
level  of  turbulence  energy  near  the  ground  surface  {z/H  <2/3).  A  possible  explanation  for 
this  underestimation  in  the  TKE  level  within  the  impingement  zone  is  that  the  linear  eddy- 
viscosity  model  used  with  the  k-e  turbulence  model  cannot  capture  the  particularly  high- 
level  of  normal-stress  anisotropy  expected  to  occur  in  this  region.  This  failure  to  properly 
capture  the  normal-stress  anisotropy  in  the  impingement  zone  will  have  two  effects  on  the 
production  of  TKE.  Firstly,  the  normal  stress  anisotropy  is  responsible  for  the  production 
of  shear  stress  — (u'rt;')  which  in  the  streamline  coordinate  system  (i.e.,  the  streamlines  and 
their  orthogonal  trajectories  define  the  orthogonal  coordinate  net)  takes  the  form  (2(u'^)  — 
{w'’^))U/R  {R  is  the  local  radius  of  streamline  curvature),  which  in  the  impingement  zone 
would  correspond  to  a  gain  in  —  (u'tt;')  owing  to  the  positive  (unstable)  flow  curvature  {R  > 
0)  in  this  region  (noting  that  2(u'^)  ^  Hence,  a  reduction  in  normal-stress  anisotropy 

would  lead  to  an  underestimation  in  —  (u'tt;'),  which  in  turn  would  lead  to  a  reduced  shear 
stress  production,  —{u'w')dU/dz,  of  the  TKE.  Secondly,  an  underestimation  of  the  normal- 
stress  anisotropy  will  lead  to  an  under-prediction  of  the  normal  stress  production  of  fc, 
—  ({u'^)  —  {w''^))dU /dx^  in  the  impingement  zone  in  front  of  the  first  building.  Note  that  in 
this  region  dU fdx  <  0  owing  to  the  deceleration  of  the  flow  in  front  of  the  windward  face 
of  the  first  building  (induced  by  the  adverse  horizontal  pressure  gradient  here),  so  the  term 
(in  streamline  coordinate  system)  — ((u'^)  —  {w^'^))dU / dx^  corresponds  to  a  gain  in  TKE. 

Figure  20  presents  a  comparison  of  measured  and  modeled  profiles  of  k  for  five  x-Iocations 
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(0.1  <  x/H  <  0.9)  over  the  first  building  rooftop.  The  model  predicts  a  monotonic  in¬ 
crease  of  the  turbulence  level  (and,  hence,  a  steady  production  of  turbulence  energy)  in  the 
streamwise  direction  along  and  near  the  rooftop  surface.  This  feature  largely  conforms  with 
the  wind  tunnel  simulations.  However,  the  maximum  TKE  level  is  underestimated  in  the 
numerical  simulations.  The  impingement  process  occuring  in  front  of  the  first  building  wall 
also  provokes  severe  concave  (unstable)  flow  curvature  in  the  region  of  the  flow  immediately 
upstream  of  the  first  building  rooftop.  It  is  expected  that  the  turbulence  kinetic  energy 
in  this  region  will  increase  significantly  due  to  the  cumulative  effect  of  concave  streamline 
curvature  as  fluid  elements  are  forced  up  and  over  the  first  building.  Indeed,  the  curvature 
term  (in  streamline  coordinates)  —{u'w')U/R  is  a  source  of  production  of  k  in  this  region 
of  positive  (unstable)  R.  Furthermore,  the  region  of  concave  streamline  curvature  extends 
over  the  first  building  rooftop.  The  subtle  effect  of  curvature  strain  on  turbulence  does  not 
appear  to  be  captured  in  the  present  model  and,  hence,  may  explain  differences  in  the  TKE 
levels  over  the  first  building  rooftop  between  the  model  predictions  and  measurements.  Al¬ 
ternatively,  the  observed  large  peak  in  the  TKE  just  above  the  first  building  rooftop  can 
be  also  interpreted  as  arising  from  an  oscillation  of  the  high  shear  layer  by  the  larger-scale 
upstream  turbulence  which  would  enhance  (and,  hence  k)  as  the  probe  at  this  level 
records  high  and  low  velocities  as  the  shear  zone  flaps  upwards  and  downwards.  The  fact 
that  such  a  large  peak  is  missing  in  the  predictions  may  suggest  that  this  process  is  not 
well  represented  in  the  k-e  turbulence  model.  Since  TKE  is  exported  downstream  by  local 
advection  and  turbulent  transport,  an  underestimation  of  the  TKE  level  in  the  region  just 
upstream  of  the  first  building  will  lead  generally  to  an  underestimation  of  the  TKE  level  all 
positions  downstream  of  the  location  of  strong  turbulence  generation. 

Figure  21  compares  computed  and  measured  vertical  profiles  of  k  within  the  first  street 
canyon  for  x-stations  in  the  range  1.1  <  x/H  <  1.9.  The  position  of  the  prominent  “nose” 
in  the  k.  profiles  lying  just  above  the  street  canyon  top  is  largely  reproduced  by  the  numerical 
simulation,  albeit  the  peak  value  of  k  is  significantly  underestimated.  This  defect  is  a 
consequence  of  the  under-prediction  of  the  TKE  level  upstream  of  this  region  owing  to  the 
inability  of  the  present  model  to  account  properly  either  for  the  effects  of  secondary  strain 
on  the  turbulence  and/or  for  the  effects  of  the  large-scale  flapping  of  the  strong  shear  layer 
at  the  canopy  top.  Just  downstream  of  the  first  building  rooftop  and  on  the  upstream 
side  of  the  first  street  canyon  (e.g.,  at  x  =  165  mm),  there  is  a  thin  intense  shear  layer 
where  the  turbulence  energy  is  large.  Moving  downstream  of  this  position,  the  evolution 
of  the  TKE  profiles  is  dominated  by  the  vertical  spreading  of  this  shear  or  mixing  layer. 
At  the  center  of  this  layer  lying  just  above  the  building  height  H,  the  peak  value  of  TKE 
attenuates  downstream,  whereas  regions  above  and  below  the  peak  “import”  turbulence 
energy  by  outward  pressure  or  turbulent  transport  from  the  core  of  the  shear  layer,  causing 
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an  increase  in  the  turbulence  energy  in  these  regions.  These  characteristics  of  the  flow  are 
largely  reproduced  by  the  numerical  simulation.  Furthermore,  the  downstream  advection 
and  vertical  diffusion  of  TKE  results  in  a  larger  level  of  turbulence  energy  being  found 
on  the  downstream  side  of  a  street  canyon  than  the  upstream  side,  a  feature  that  is  also 
present  in  the  model  predictions.  Within  the  street  canyon,  minimum  levels  of  TKE  are 
found  on  the  upstream  side  of  the  canyon  (protected  region)  where  the  TKE  is  essentially 
independent  of  height.  These  features  are  reproduced  in  the  simulations,  although  the  near 
constant  TKE  level  on  the  upstream  side  of  the  canyon  is  underestimated  by  the  model. 
Finally,  it  appears  that  the  rate  of  vertical  spreading  of  the  shear  or  mixing  layer  in  the 
numerical  simulation  is  smaller  than  in  the  measurements.  This  discrepancy  may  be  the 
result  of  an  underestimation  of  the  turbulent  viscosity  Ut  [cf-  Equation  (5)]  owing  to  the 
smaller  predicted  TKE  levels. 

Figures  22  and  23  exhibit  the  comparison  of  the  predicted  values  of  TKE  with  experimental 
data,  respectively,  over  the  second  building  rooftop  and  the  second  street  canyon  for  the 
same  x-stations  as  the  mean  velocity  field  described  above.  In  contrast  to  the  behavior 
over  the  first  building  rooftop,  the  peak  turbulence  energy  decreases  monotonically  (ap¬ 
proximately  or  better)  as  one  moves  downstream  across  the  second  building  rooftop.  Such 
a  feature  was  not  found  in  the  numerical  simulations.  This  observed  property  of  the  flow 
may  be  the  result  of  streamline  curvature  on  TKE  production.  Specifically,  the  reduced 
production  (loss)  of  TKE  along  the  second  building  rooftop  arises  from  two  effects:  (1) 
the  TKE  production  term  due  to  curvature  strain,  -{u'w^)U /R,  represents  a  (negative) 
loss  term  in  this  region  of  negative  (stable)  R  (viz.,  the  downward  or  convex  curvature 
of  the  mean  streamlines  over  the  second  rooftop  suppresses  the  turbulence  here);  and,  (2) 
the  reduction  of  the  shear  stress  -{u'w')  due  to  convex  streamline  curvature  will  lead  to 
reduced  shear  production.  Again,  it  is  seen  that  the  subtle  effect  of  curvature  strain  on  the 
turbulence  is  not  adequately  captured  by  the  present  model.  All  the  qualitative  features 
of  the  turbulence  energy  profiles  are  predicted  moderately  well  in  the  third  street  canyon, 
although  the  TKE  levels  are  underestimated  at  all  vertical  levels.  However,  the  degree  of 
underestimation  is  not  as  great  as  it  was  for  the  first  street  canyon. 

Comparisons  of  computed  and  measured  profiles  of  k  over  the  third,  fourth,  fifth,  and  sixth 
building  rooftops  and  street  canyons  are  presented  in  Figures  24  to  31.  Figure  32  compares 
computed  TKE  profiles  with  measured  ones  over  the  seventh  building  rooftop.  Unlike 
the  mean  velocity  field  which  achieves  streamwise  equilibrium  in  the  array  interior  when 
the  flow  reaches  the  third  building,  the  TKE  never  seems  to  reach  streamwise  equilibrium 
and,  consequently,  exhibits  downstream  development  throughout  the  entire  array  interior. 
Interestingly,  the  level  of  TKE  over  the  rooftops  in  the  array  interior  change  negligibly  as 
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one  moves  in  the  streamwise  direction  across  the  rooftops,  suggesting  that  there  is  a  local 
equilibrium  established  here  between  the  production  and  dissipation  of  TKE.  In  contrast, 
the  peak  TKE  levels  decrease  monotonically  as  one  moves  downstream  from  the  upstream 
to  the  downstream  side  of  the  street  canyons  as  the  TKE  in  the  shear  layer  above  the 
building  top  is  “exported”  by  turbulent  diffusion  into  the  canyons.  Furthermore,  the  TKE 
levels  in  street  canyons  further  downstream  are  successively  smaller  owing  to  the  fact  that 
the  peak  TKE  in  the  mixing  layer  above  each  canyon  attenuates  downstream  (and,  hence, 
there  is  less  TKE  available  to  transport  into  each  successive  street  canyon  downstream). 
These  qualitative  features  of  the  flow  are  largely  reproduced  by  the  numerical  simulation. 
The  quantitative  agreement  between  the  computed  and  measured  TKE  profiles  generally 
improve  with  increasing  downstream  distance. 

Figures  33  and  34  compare  computed  and  measured  TKE  profiles  in  the  exit  region  of 
the  array  directly  downstream  of  the  last  building  covering  the  range  13.1  <  x/F  <  20.5. 
Overall,  the  TKE  profiles  in  the  exit  region  appear  to  be  very  similar  in  shape,  and  all 
the  qualitative  features  in  these  profiles  appear  to  be  correctly  predicted;  the  quantitative 
agreement  is  also  fairly  good.  Just  downstream  of  the  last  obstacle  and  below  the  build¬ 
ing  height  at  X  =  1965  and  1988  mm  (cf.  Figure  33),  the  turbulence  energy  is  practically 
eliminated  within  this  sheltered  “quiet”  zone,  a  feature  that  is  correctly  reproduced  by  the 
model.  For  the  profiles  in  the  range  from  2025  to  2775  mm,  the  position  of  the  maximum 
TKE  is  correctly  predicted.  This  position  corresponds  to  the  center  of  the  curved  shear  layer 
that  bounds  or  envelops  the  recirculation  zone  behind  the  last  building.  The  characteris¬ 
tics  of  the  vertically  spreading  mixing  zone  (or,  wake  zone)  above  and  downstream  of  the 
recirculation  zone  are  predicted  well  by  the  model.  Finally,  the  rate  of  recovery  of  the  TKE 
in  the  far  wake  region  to  the  upstream  reference  condition  is  reproduced  correctly.  Even  at 
X  =  3075  mm  or,  equivalently,  x/F  =  7.5  downstream  of  the  last  building  (measurement 
furthest  downstream),  the  TKE  profile  has  the  appearance  of  a  detached  mixing  layer  de¬ 
spite  the  fact  that  the  mean  flow  has  re-attached  at  about  x/F  =  3.5  downstream  of  the 
last  building.  This  profile  form  is  a  residue  of  the  strong  shear  or  mixing  layer  that  detaches 
from  the  last  building  rooftop  and  spreads  outward  by  pressure  or  turbulent  diffusion  until 
this  vertically  spreading  mixing  layer  makes  contact  with  the  ground  surface.  This  is  an 
example  of  disequilibrium  that  results  in  the  TKE  profile  at  even  the  furthest  downstream 
measurement  position  not  having  recovered  to  the  reference  far  upstream  state. 
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Summary  cmd  conclusions 


A  linear  eddy-viscosity  k-e  turbulence-closure  model,  applied  with  wall  functions,  has  been 
used  to  calculate  the  disturbed  mean  flow  and  turbulence  near  and  about  an  array  of  two- 
dimensional  buildings.  The  primary  purpose  of  this  exercise  was  to  assess  critically  the 
model’s  predictive  performance  by  comparison  of  the  mean  flow  and  turbulence  statistics 
with  a  very  comprehensive  and  detailed  wind  tunnel  simulation.  It  should  be  emphasized 
that  for  this  numerical  simulation,  the  closure  constants  in  the  k-e  turbulence  model  re¬ 
tained  their  original  (standard)  values.  In  consequence,  all  model  predictions  were  obtained 
without  adjustment  any  of  the  closure  constants.  This  ensures  that  the  process  did  not 
simply  degenerate  into  a  “curve-fitting”  exercise  . 

The  comparison  of  predicted  and  measured  streamwise  and  vertical  velocity  profiles  indi¬ 
cated  that  most  of  the  qualitative  features  of  the  mean  flow  field — such  as  speed-up  over 
the  first  building  rooftop,  the  presence  of  a  thin  separation  layer  near  the  surface  of  the  first 
rooftop  and  the  absence  of  such  a  feature  in  all  the  succeeding  rooftops,  the  standing  vortex 
in  the  street  canyon,  the  upward  vertical  motion  on  the  first  rooftop,  the  downward  vertical 
motion  on  the  second  rooftop,  the  thin  strong  shear  layer  that  forms  at  building  height, 
the  recirculation  zone  behind  the  last  building,  and  the  development  of  the  mixing  layer 
above  and  downstream  of  this  recirculation  zone,  and  the  rate  of  momentum  recovery  in 
the  far  wake  region — are  correctly  predicted.  The  quantitative  agreement  is  also  fairly  good 
throughout  most  of  the  flow  domain,  with  the  greatest  discrepancy  in  the  mean  velocity  field 
arising  from  the  prediction  of  a  smaller  than  observed  separation  zone  in  the  impingement 
region  upstream  of  the  first  building  and  an  underestimation  of  the  magnitude  of  both  the 
downward  vertical  velocity  on  the  downstream  side  and  the  upward  vertical  velocity  on  the 
upstream  side  of  a  street  canyon. 

The  turbulence  kinetic  energy  from  the  numerical  and  physical  simulations  show  striking 
resemblances  and,  indeed,  most  of  the  qualitative  features  of  the  turbulence  energy  as 
its  develops  near  and  through  the  array  of  buildings  is  adequately  predicted.  In  terms 
of  quantitative  agreement,  it  is  found  that  the  TKE  level  is  generally  underestimated  by 
the  model.  This  underestimation  may  be  caused  by  the  inability  of  the  model  to  capture 
the  disproportionate  and  inhomogeneous  sensitivity  of  the  Reynolds  stress  components  to 
secondary  strains  in  the  flow,  particularly  those  associated  with  curvature.  Specifically, 
the  production  of  TKE  just  upstream  and  above  the  first  building  and  across  the  first 
building  rooftop  arising  from  the  mean  streamline  concave  (destabilising)  curvature  strain 
acting  through  the  shear  stress  is  not  adequately  captured  by  the  model.  Because  the 
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TKE  is  exported  downstream  and  vertically  by  advection  and  by  turbulent  and /or  pressure 
diffusion,  the  under-prediction  of  the  TKE  at  and  near  the  leading  edge  of  the  first  building 
(i.e.,  the  location  where  the  turbulence  energy  production  is  greatest)  will  result  in  reduced 
TKE  levels  at  all  positions  downstream  of  this  location. 

The  k-e  turbulence-closure  model  with  an  isotropic,  linear  eddy  viscosity  is  the  simplest 
complete  turbulence  model  (viz.,  no  advance  knowledge  of  any  property  of  the  turbulence 
is  required  for  the  simulation  other  than  the  initial  and/or  boundary  conditions  for  the 
problem)  that  is  available  currently,  and  its  moderately  good  predictive  performance  of 
the  complicated  developing  flow  through  an  array  of  two-dimensional  buildings  without  the 
need  to  adjust  any  closure  constants  is  encouraging.  This  model  may  be  useful  as  a  general- 
purpose  simulator  of  urban  flows  since  it  is  simple  enough  to  be  tractable  numerically  and, 
hence,  not  require  excessive  computing  time.  Notwithstanding  certain  weaknesses  of  the  k-e 
turbulence  model  for  urban  flows  as  discussed  in  this  report,  the  model  remains  attractive 
because  of  its  formalistic  simplicity,  numerical  robustness,  and  computational  economy. 
How'ever,  before  this  model  (or,  variants  of  it)  can  be  used  for  this  purpose,  it  needs  to 
be  validated  against  more  experimental  data  for  flow  in  various  obstacle  arrays.  From  this 
perspective,  the  availability  of  extensive  and  accurate  experimental  data  for  flow  through 
various  obstacle  arrays  is  urgently  required. 
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Figure  1.  The  two-dimensional  building  array  used  in  the  wind  tunnei  experiment 
The  array  consisted  of  rectanguiar  biocks  with  an  equal  height  and  length  of  150 
mm.  with  each  block  extending  waii-to-wall  across  the  wind  tunnel  in  the 
spanwise  direction. 
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Figure  2.  Computational  grid  used  for  numerical  simulation  of  the  two- 
dimensional  building  array. 
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Figure  3.  Vertical  profiles  of  the  mean  streamwise  velocity  U  at  specified  x- 
locations  in  the  impingement  zone  (-2.5  <=  x/H  <=  -0.1)  upstream  of  the  front 
face  of  the  first  building  compared  with  time-averaged  wind  tunnel 
measurements. 
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Figure  4.  Vertical  profiles  of  the  mean  streamwise  velocity  U  at  specified  x- 
locations  over  the  roof  (0.1  <=  xlH  <=  0.9)  of  the  first  building  compared  with 
time-averaged  wind  tunnel  measurements. 
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Figure  5.  Vertical  profiles  of  the  mean  streamwise  velocity  U  at  specified  x- 
locations  within  the  first  urban  canyon  (1 .1  <=  x/H  <=  1 .9)  compared  with  time- 
averaged  wind  tunnel  measurements. 
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Figure  6.  Vertical  profiles  of  the  mean  streamwise  velocity  U  at  specified  x- 
locations  over  the  roof  (2.1  <=  x/H  <=  2.9)  of  the  second  building  compared 
with  time-averaged  wind  tunnel  measurements. 
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Figure  7.  Vertical  profiles  of  the  mean  streamwise  velocity  U  at  specified  x- 
locations  within  the  second  urban  canyon  (3.1  <=  x/H  <=  3.9)  compared  with 
time-averaged  wind  tunnel  measurements. 
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Figure  8.  Vertical  profiles  of  the  mean  streamwise  velocity  U  at  specified  x- 
locations  over  the  roof  (4.1  <=  x/H  <=  4.9)  of  the  third  building  compared  with 
time-averaged  wind  tunnel  measurements. 
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Figure  9.  Vertical  profiles  of  the  mean  streamwise  velocity  U  at  specified  x- 
locations  in  the  exit  region  (13.1  <=  x/H  <=  14.5)  downstream  of  the  back  face 
of  the  last  building  compared  with  time-averaged  wind  tunnel  measurements. 
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Figure  10.  Vertical  profiles  of  the  mean  streamwise  velocity  U  at  specified  x- 
locations  in  the  exit  region  (15.5  <=  x/H  <=  20.5)  downstream  of  the  back  face 
of  the  last  building  compared  with  time-averaged  wind  tunnel  measurements. 
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Figure  11.  Vertical  profiles  of  the  mean  vertical  velocity  W  at  specified  x-locations 
in  the  impingement  zone  (-2.5  <=  x/H  <=  -0.1)  upstream  of  the  front  face  of  the 
first  building  compared  with  time-averaged  wind  tunnel  measurements. 
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Figure  12.  Vertical  profiles  of  the  mean  vertical  velocity  Wat  specified x-locations 
over  the  roof  (0.1  <=  x/H<=  0.9)  of  the  first  building  compared  with  time- 
averaged  wind  tunnel  measurements. 
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Figure  13.  Vertical  profiles  of  the  mean  vertical  velocity  W  at  specified  x-locations 
within  the  first  urban  canyon  (1.1  <=  x/H  <=  1 .9)  compared  with  time-averaged 
wind  tunnel  measurements. 
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Figure  14.  Vertical  profiles  of  the  mean  vertical  velocity  Wat  specified  x-locations 
over  the  roof  (2.1  <=  x/H  <=  2.9)  of  the  second  building  compared  with  time- 
averaged  wind  tunnel  measurements. 
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Figure  15.  Vertical  profiles  of  the  mean  vertical  velocity  \N  at  specified  x-locations 
within  the  second  urban  canyon  (3.1  <=  x/H  <=  3.9)  compared  with  time- 
averaged  wind  tunnel  measurements. 
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Figure  16.  Vertical  profiles  of  the  mean  vertical  velocity  W  at  specified  x-locations 
over  the  roof  (4.1  <=  x/H  <=  4.9)  of  the  third  building  compared  with  time- 
averaged  wind  tunnel  measurements. 


36 


DRESTR  2001-1 15 


W<m/s) 


Figure  17.  Vertical  profiles  of  the  mean  vertical  velocity  W  at  specified  x-locations 
in  the  exit  region  (13.1  <=  x/H  <=  14.5)  downstream  of  the  back  face  of  the  last 
building  compared  with  time-averaged  wind  tunnel  measurements. 
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Figure  18.  Vertical  profiles  of  the  mean  vertical  velocity  W  at  specified  x-locations 
in  the  exit  region  {1 5.5  <=  x/H  <=  20.5)  downstream  of  the  back  face  of  the  last 
building  compared  with  time-averaged  wind  tunnel  measurements. 


Figure  19.  Vertical  profiles  of  the  turbulence  kinetic  energy  k  at  specified  x- 
locations  in  the  impingement  zone  (-2.5  <=  x/H  <=  -0.1)  upstream  of  the  front 
face  of  the  first  building  compared  with  time-averaged  wind  tunnel 
measurements. 
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Figure  20.  Vertical  profiles  of  the  turbulence  kinetic  energy  k  at  specified  x- 
locations  over  the  roof  (0.1  <=  x/H  <=  0.9)  of  the  first  building  compared  with 
time-averaged  wind  tunnel  measurements. 
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Figure  21.  Vertical  profiles  of  the  turbulence  kinetic  energy  k  at  specified  x- 
locations  within  the  first  urban  canyon  (1.1  <=  xlH  <=  1 .9)  compared  with  time- 
averaged  wind  tunnel  measurements. 
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Figure  22.  Vertical  profiles  of  the  turbulence  kinetic  energy  at  specified  x- 
locations  over  the  roof  (2.1  <=  x/H  <=  2.9)  of  the  second  building  compared 
with  time-averaged  wind  tunnel  measurements. 
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Figure  23.  Vertical  profiles  of  the  turbulence  kinetic  energy  k  at  specified  x- 
locations  within  the  second  urban  canyon  (3.1  <=  x/H  <=  3.9)  compared  with 
time-averaged  wind  tunnel  measurements. 
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Figure  24.  Vertical  profiles  of  the  turbulence  kinetic  energy  k  at  specified  x- 
locations  over  the  roof  (4.1  <=  a//-/  <=  4.9)  of  the  third  building  compared  with 
time-averaged  wind  tunnel  measurements. 


44 


DRESTR  2001-1 15 


Figure  25.  Vertical  profiles  of  the  turbulence  kinetic  energy  k  at  specified  x- 
locations  within  the  third  urban  canyon  (5.1  <=  x/H  <=  5.9)  compared  with  time 
averaged  wind  tunnel  measurements 
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Figure  26.  Vertical  profiles  of  the  turbulence  kinetic  energy  k  at  specified  x- 
locations  over  the  roof  (6.1  <=  x/H  <=  6.9)  of  the  fourth  building  compared  with 
time-averaged  wind  tunnel  measurements. 
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Figure  27.  Vertical  profiles  of  the  turbulence  kinetic  energy  k  at  specified  x- 
locations  within  the  fourth  urban  canyon  (7.1  <=  x/H  <=  7.9)  compared  with 
time-averaged  wind  tunnel  measurements. 
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Figure  28.  Vertical  profiles  of  the  turbulence  kinetic  energy  k  at  specified  x- 
locations  over  the  roof  (8.1  <=  x/H  <=  8.9)  of  the  fifth  building  compared  with 
time-averaged  wind  tunnel  measurements. 
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Figure  29.  Vertical  profiles  of  the  turbulence  kinetic  energy  k  at  specified  x- 
locations  within  the  fifth  urban  canyon  (9.1  <=  x/H  <=  9.9)  compared  with  time- 
averaged  wind  tunnel  measurements. 
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Figure  30.  Vertical  profiles  of  the  turbulence  kinetic  energy  k  at  specified  x- 
locations  over  the  roof  (10.1  <=  x/H  <=  10.9)  of  the  sixth  building  compared 
with  time-averaged  wind  tunnel  measurements. 
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Figure  31.  Vertical  profiles  of  the  turbulence  kinetic  energy  k  at  specified  x- 
locations  within  the  sixth  urban  canyon  (11.1  <=  x/H  <=11 ,9)  compared  with 
time-averaged  wind  tunnel  measurements. 
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Figure  32.  Vertical  profiles  of  the  turbulence  kinetic  energy  k  at  specified  x- 
locations  over  the  roof  (12.1  <=  x/H  <=  ^2.9)  of  the  seventh  building  compared 
with  time-averaged  wind  tunnel  measurements. 
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Figure  33.  Vertical  profiles  of  the  turbulence  kinetic  energy  k  at  specified  x- 
locations  in  the  exit  region  (13.1  <=  x/H  <=  14.5)  downstream  of  the  back  face 
of  the  last  building  compared  with  time-averaged  wind  tunnel  measurements. 
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Figure  34.  Vertical  profiles  of  the  turbulence  kinetic  energy  k  at  specified  x- 
locations  in  the  exit  region  (15.5  <=  x/H  <=  20.5)  downstream  of  the  back  face 
of  the  last  building  compared  with  time-averaged  wind  tunnel  measurements. 
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